Plotting a tange

library(plotly)
library(pracma)
library(Deriv)
plot_tangent_surface = function(f, x_0, y_0, grid_x = -10, grid_y = 10) {
  # Create meshgrid for plotting 
  m = meshgrid(grid_x:grid_y)
  nrows = nrow(m$X)
  X = as.vector(m$X)
  Y = as.vector(m$Y)
  z_0 = f(x_0, y_0)
  df_dx = Deriv(f, 'x')
  df_dy = Deriv(f, 'y')
  
  F_MAT = f(X, Y)
  SURFACE = z_0 + df_dx(x_0, y_0)*(X - x_0) + df_dy(x_0, y_0)*(Y - y_0)
  
  F_MAT = matrix(F_MAT, nrow = nrows)
  SURFACE = matrix(SURFACE, nrow = nrows)
  
  p <- plot_ly(z = ~F_MAT) %>% add_surface(
    contours = list(
        z = list(
        show=TRUE,
        usecolormap=TRUE,
        highlightcolor="#ff0000",
        project=list(z=TRUE)
      )
    )
  ) %>%
    layout(
      scene = list(
        camera=list(
          eye = list(x=1.87, y=0.88, z=-0.64)
        )
      )
    ) %>%
    add_surface(z = ~SURFACE)
  
  p
}
# Define f
f = function(x, y) {
  x^2  + y^2
}
x_0 = -3
y_0 = -3
# debug(plot_tangent_surface)
plot_tangent_surface(f, x_0, y_0, -15, 15)
LS0tCnRpdGxlOiAiVmVjdG9yIENhbGN1bHVzOiBQbG90dGluZyBUYW5nZW50IFN1cmZhY2UiCm91dHB1dDogaHRtbF9ub3RlYm9vawplZGl0b3Jfb3B0aW9uczogCiAgY2h1bmtfb3V0cHV0X3R5cGU6IGNvbnNvbGUKLS0tCgpQbG90dGluZyBhIHRhbmdlCgpgYGB7cn0KbGlicmFyeShwbG90bHkpCmxpYnJhcnkocHJhY21hKQpsaWJyYXJ5KERlcml2KQoKcGxvdF90YW5nZW50X3N1cmZhY2UgPSBmdW5jdGlvbihmLCB4XzAsIHlfMCwgZ3JpZF94ID0gLTEwLCBncmlkX3kgPSAxMCkgewogICMgQ3JlYXRlIG1lc2hncmlkIGZvciBwbG90dGluZyAKICBtID0gbWVzaGdyaWQoZ3JpZF94OmdyaWRfeSkKICBucm93cyA9IG5yb3cobSRYKQogIFggPSBhcy52ZWN0b3IobSRYKQogIFkgPSBhcy52ZWN0b3IobSRZKQoKCiAgel8wID0gZih4XzAsIHlfMCkKCiAgZGZfZHggPSBEZXJpdihmLCAneCcpCiAgZGZfZHkgPSBEZXJpdihmLCAneScpCiAgCiAgRl9NQVQgPSBmKFgsIFkpCiAgU1VSRkFDRSA9IHpfMCArIGRmX2R4KHhfMCwgeV8wKSooWCAtIHhfMCkgKyBkZl9keSh4XzAsIHlfMCkqKFkgLSB5XzApCiAgCiAgRl9NQVQgPSBtYXRyaXgoRl9NQVQsIG5yb3cgPSBucm93cykKICBTVVJGQUNFID0gbWF0cml4KFNVUkZBQ0UsIG5yb3cgPSBucm93cykKICAKICBwIDwtIHBsb3RfbHkoeiA9IH5GX01BVCkgJT4lIGFkZF9zdXJmYWNlKAogICAgY29udG91cnMgPSBsaXN0KAogICAgICAgIHogPSBsaXN0KAogICAgICAgIHNob3c9VFJVRSwKICAgICAgICB1c2Vjb2xvcm1hcD1UUlVFLAogICAgICAgIGhpZ2hsaWdodGNvbG9yPSIjZmYwMDAwIiwKICAgICAgICBwcm9qZWN0PWxpc3Qoej1UUlVFKQogICAgICApCiAgICApCiAgKSAlPiUKICAgIGxheW91dCgKICAgICAgc2NlbmUgPSBsaXN0KAogICAgICAgIGNhbWVyYT1saXN0KAogICAgICAgICAgZXllID0gbGlzdCh4PTEuODcsIHk9MC44OCwgej0tMC42NCkKICAgICAgICApCiAgICAgICkKICAgICkgJT4lCiAgICBhZGRfc3VyZmFjZSh6ID0gflNVUkZBQ0UpCiAgCiAgcAp9CgojIERlZmluZSBmCmYgPSBmdW5jdGlvbih4LCB5KSB7CiAgeF4yICArIHleMgp9Cgp4XzAgPSAtMwp5XzAgPSAtMwoKIyBkZWJ1ZyhwbG90X3RhbmdlbnRfc3VyZmFjZSkKcGxvdF90YW5nZW50X3N1cmZhY2UoZiwgeF8wLCB5XzAsIC0xNSwgMTUpCmBgYAo=